Astronomy & Astrophysics manuscript no. 161 14 


©ESO 2013 


January 1 1, 2013 





o 

(N 



o 



(N 
> 

\o 
l> 
o 

o 



X 



Rotating massive main-sequence stars II: 
Simulating a population of LMC early B-type stars as a test of 

rotational mixing 

I. C. J. Evans^ I. Hunter^, A. de Kotei-^ ', N. Langer^ ', P. L. Dufton^, M. Cantiello^ C. Trundle'*, 

D. J. Lennon^ S.E. de Mink^-^ S.-C. Yoon^, and P. Anders^ 

' Astronomical Institute, Utrecht University, Princetonplein 5, 3584 CC, Utreclit, The Netherlands 

- University of Vienna, Department of Astronomy, Tiirkenschanzstr. 17, A-1180 Vienna, Austria 

^ UK Astronomy Technology Centre, Royal Observatory Edinburgh, Blackford Hill, Edinburgh, EH9 3HJ, UK 

Astrophysics Research Centre, School of Mathematics & Physics, The Queens University of Belfast, Belfast, BT7 INN, Northern 

Ireland, UK 

' Astronomical Institute Anton Pannekoek, University of Amsterdam, Kruislaan 403, 1098 SJ, Amsterdam, The Netherlands 
* Argelander-Institut fiir Astronomic der Universitat Bonn, Auf dem Hiigel 71, 53121 Bonn, Germany 
' Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA 
** Hubble Fellow 



Received: / Accepted: 

ABSTRACT 

Context. Rotational mixing in massive stars is a widely applied concept, with far-reaching consequences for stellar evolution, nucle- 
osynthesis, and stellar explosions. 

Aims. Nitrogen surface abundances for a large and homogeneous sample of massive B-type stars in the Large Magellanic Cloud 
(LMC) have recently been obtained by the ESO VLT-FLAMES Survey of Massive Stars. This sample is the first to cover a broad 
range of projected stellar rotational velocities, with a large enough sample of high quality data to allow for a statistically significant 
analysis. Here we use the sample to provide the first rigorous and quantitative test of the theory of rotational mixing in massive stars. 
Methods. We calculated a grid of stellar evolution models, using the VLT-FLAMES sample to calibrate some of the uncertain mixing 
processes. We developed a new population-synthesis code, which uses this grid to simulate a large population of stars with masses, 
ages, and rotational velocity distributions consistent with those from the VLT-FLAMES sample. The synthesized population is then 
filtered by the selection effects in the observed sample, to enable a direct comparison between the empirical results and theoretical 
predictions. 

Results. Our simulations reproduce the fraction of stars without significant nitrogen enrichment. However, the predicted number of 
rapid rotators with enhanced nitrogen is about twice as large as found observationally. Furthermore, two groups of stars (one consist- 
ing of slowly rotating, nitrogen-enriched objects and another consisting of rapidly rotating un-enriched objects) cannot be reproduced 
by our single-star population synthesis. 

Conclusions. Physical processes in addition to rotational mixing appear to be required to understand the population of massive main- 
sequence stars from the VLT-FLAMES sample. We discuss the possible role of binary stars and magnetic fields in the interpretation 
of our results. We find that the population of slowly rotating nitrogen-enriched stars is unlikely to be produced via mass transfer and 
subsequent tidal spin-down in close binary systems. A conclusive assessment of the role of rotational mixing in massive stars requires 
a quantitative analysis that also accounts for the effects of binarity and magnetic fields. 

Key words. Astronomical methods: statistical - Stars: abundances, evolution, early-type, rotation, massive 

1 . Introduction tators viewed nearly pole-on. Since the range of observed enrich- 

ments was consistent with evolutionary predictions of rapidly 
Observations of chemi c ally enriched rn a in-sequence stars (e.g . rotating stars, this wa s considered as confirrning rotational mix- 
Gies & Lambert '1992!; iHerrerol Il993l IVrancken et aP l2000t ing in massive stars dHeger & Lan gel^ l2000^ iMevnet & Maedel^ 



Przvbilla et al..,2010i) have fostered the idea that rotationally in- I2000h . 



duced mixing may bring fusion products from the core of mas- It was only in the fr amework of the E SO VLT-FLAMES 
sive stars to their surface during core hydrogen burning. While Survey of Massive Stars (lEvans et al.ll2005h that nitrogen abun- 
the potential consequences of rotational mixing fo r the evolution dances were determined for massive main- sequence stars with 
of massive stars have been shown to be dramatic (lMaedeill2000t a wide range of p rojected rotation rates (iHunter et al.l l2007t 
lYoon & Lange3l2005h . direct tests of the models have proven to iTrundle et al.ll2007h . Further analysis of the early B-type stars 
be difficult. A major obstacle of such tests is that observations in the Large Magellanic Cloud (LMC) from this survey re- 
have typically been restricted to stars that have low projected ro- vealed the existence of a significant population of nitrogen 
tation rates, to minimize any rotational broadening of the stellar enhanced stars with int rinsically slow rot ation (iHunter et al.l 
absorption lines. However, the often small sample size has al- I2008al) . The findings of iMorel etall (l2006b argue for a simi- 
lowed speculation that nitrogen-enriched stars might be rapid ro- lar population in our Galaxy. The origin of the nitrogen en- 
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hancement in these stars is not understood, but, since they are 
slow rotators, it seems doubtful that rotational mixing can ex- 
plain it. It is therefore conceivable that samples of enriched 
stars, with apparently slow rotation velocities (which substan- 
tiated the most convincing direct evidence for rotational mix- 
ing in massive main-sequence stars over the past decades), may 
instead be re lated to c ompletely different physical processes - 
e.g. b inarity dLanger et a l. 2008) or magnetic fields (Morel et all 
I2OO8I) . In this context, it is worth noting that a wealth of 
atmospheric and wind properties of massive main-sequence 
stars are yet to be understood, such as the winds of massive 
stars with lum inosities belo w ~10^'^Lo (e.g'Bouret et al. 2OOI 
iMartins et al]l2004 . 2005; Mokiem et al.ll2b07; Puis et al. 2Q(M, 
intermitten t discrete absorption components in UV absorp- 
tion lines dPrinia & HowarthllI988l: iKaper et al.llI997[), micro - 
and macro-turbulence (ICantiello et al.ll2009t lAerts et al.ll2009h . 
and non-ther mal X-ray emission (iBabel & Montmerld Il997t 
lud-Doula &"bwocki2002<) . 

Putting the probl em of the enriched slow r otators aside, both 
iHunteret all ( 12008 ah and lMaeder et al.l ( |2009|) suggested that the 
majority of the newly discovered population of rapidly rotating, 
nitrogen enhanced main-sequence stars in the VLT-FLAMES 
survey may agree with the predictions of rotational mixing. 
However, their results remain ambiguous, as the existence of a 
population of rapidly rotating, nitrogen enriched, apparently sin- 
gle stars is also a prime predic tion of close binary e v olution mod- 
els ac counting for rotation (iPetrovic et al.l l2005t iLanger et alj 
120081) . which is independent of the mechanism of rotational mix- 
ing. 

In this context, a quantitative effort to reproduce the physical 
properties of the LMC B-type sample from the VLT-FLAMES 
survey (see Sec.|2]i using stellar evolution models for single stars 
which allow for the effects of rotational mixing, is a logical next 
step. A future step will be to also account for binary evolution. 
To achieve the present aim, we use a dense grid of stel lar evolu- 
tion models for main-sequence stars (Brott et ani2010[ Paper I), 
for which the initial abundances, convective core overshooting, 
and rotational mixing were calibrated using the results from the 
VLT-FLAMES survey. These models are then incorporated in a 
newly developed population-synthesis code (Sec. [3]), which we 
employ in an attempt to reproduce the properties of the LMC 
sample of B-type stars from the VLT-FLAMES survey. The ob- 
served sample is briefly discussed in Sec. |2] Our results are de- 
scribed in Sec.|4] We draw our conclusions in Sec.|5] 



2. Observational sample 

To investigate the predictions of rotational mixing, in this paper 
we focus on population synthesis of the 107 main-sequence B - 
type stars in the LMC discussed bv lHunter eTaP (I2008ail2009ah . 
These were taken from the VLT-FLA MES obs ervations in the 
two LMC fields - Nl 1 and NGC 2004 dEvans et al. 2006). Nil 
is the second largest Hii region in the LMC, with multiple 
clumps of star formation in which the VLT-FLAMES obser- 
vations sample a range of stellar ages and spatial structures. 
NGC 2004 is an older, fairly condensed cluster, in which the 
targets were mostly in the outer region (the core is too densely 
populated for the FLAMES-Medusa fibres) and out into the sur- 
rounding field population. A total of 225 O- and early B-type 
stars were observed with VLT-FLAMES in the two LMC point- 
ings. 

To constrain our synthesis models we need to understand 
the selection effects which influenced the original observed sam- 
ple and the subsample considered bv lHunter et al.l (l2008al) . The 



most significant factor in the observed sample was the faint cut- 
off {V < 15.5™, to ensure sufficient signal-to-noise in the final 
spectra), combined with functional issues such as crowding of 
potential targets. 

In the NGC 2004 fi eld we note that, contrary to the statement 
bv lEvans et alJ d2p06l) . some of the Be-type stars identified by 
iKelleret al.l ((19991) ' were actually hardwired out of the input cat- 
alog used for target selection. This could potentially bias the final 
distribution of Be- relative to B-type stars in this field. However, 
given the magnitude (and color) constraints on targets, only 18 
such st ars would have been i ncluded in our catalog of potential 
targets. iHunter et al.l (l200 8b!) investigated the potential for sim- 
ilar selection effects in NGC 330 (in the SMC) using the Fibre 
Positioner Observation Support Software (FPOSS) with the in- 
clusion of previously excluded stars. We have adopted a similar 
method for NGC 2004, running FPOSS ten times to compare the 
resulting configura tions. Of the previously excluded stars from 
iKeller et aH dl999h . six stars were included in the resulting fibre 
configuration in two of the test runs, seven stars were included 
in six of the tests, and eight stars were included in the remain- 
ing two tests. Of course, by including some of these previously 
excluded stars, a corresponding number of the stars actually ob- 
served would not have been included, some of which might well 
be Be-type stars (18 of the 116 stars observed with the Medusa 
fibres in NGC 2004 are Be-type, i.e. 15%). There may be a weak 
selection effect relating to the observed sample of Be-type stars, 
but we believe that this does not unduly bias the overall sam- 
ple in NGC 2004 - other effects such as the spatial distribution 
of targets and the magnitude/color cuts would have been more 
significant. In short, the survey provided a relatively unbiased 
sample of the bright (V < 15.5™) main-sequence stars in the 
observed regions. 



2.1. Stellar parameters 

Physical properties of the LMC B-type stars were determined 
in a series of papers. The stellar parameters used in our 
LMC popula t ion syn thesis (Nil and NGC 2004) are those from 
iHunter et al.l d2008bl) . which incorporated r esult s for narrow- 
lined B-type stars from IHunter et al.l (l2007[) and Trundle et al.l 
d2007l) . (and for O-type stars from lMokiem et al.ll2006l) . 

Physical parameters and abundances for the B-type stars 
were obtai ned using the grid o f non-LTE model atmospheres de- 
scribed by Dufton et al. (2005'), calculated using the tlusty code 
(iHubenv & Lanz 1995) . Effective temperatures were estimated 
from the silicon or helium spectrum, while in those stars where 
adequate spectral diagnostics were not available, temperatures 



were adopted on the b asis of their spectral classifications us- 
ing the calibration from iTrundle et aL l (l2007l) . Surface gravities 
{\ogg) were deduced from the hydrogen line profiles. Projected 
rotational velocit ies (usin/) were det ermined using the profile- 
fitting methods o f iDufton etal]d2006l) , i.e. using the helium lines 
(primarily He I at /14026 A) in the majority of stars and using 
metal lines (such as Mgii, Sim) at low ysin/. Surface abun- 
dance s of C, N, O, Mg, and Si were also de termined where pos- 
sible (ITrundle et al.ll2007l: IHunter et al.ll200 7. 2008a b, 2009a). 

Results for 107 main-sequence B-type stars were given by 
IHunter et al.l d2009al) . From the original sample of 225 O- and 
early B-type stars: 41 were O-type; 1 was a W-R star; 24 were 
Be-type (for which N abundances were not determined); 23 ad- 
ditional B-type stars were not analyzed due to doubled-lined 
binarity, 'shell-like' emission, and low-quality spectra; 29 stars 
have log § < 3.2. Note however that the observations of the VLT- 
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Flames survey were not optimized for detection of binaries, so 
some almost certainly remain undetected. In Appendix |A] we 
provide mass and velocity distributions of the observational sam- 
ple. We refer to Sec. l3.3l for a description of how the selection 
effects are applied to the population synthesis model. 



3. Stellar models and population synthesis 

We have developed the population synthesis code 
STARMAKER to study the statistical properties of early- 
type stars for a specified initial rotational velocity, initial mass 
function (IMF), and given star formation history (SFH). This 
method can be applied to quantitatively constrain the effects 
of rotation (and mass loss) on stellar evolution if the dataset to 
which its results are compared is sufficiently large to provide 
good sampling of the (random) orientation of rotation axes, and 
of the velocity and age distribution. On a given grid of stellar 
evolution models, STARMAKER interpolates model stars of a 
specified initial mass, age and initial rotational velocity within 
the grid. In addition to the main stellar parameters (effective 
temperature, luminosity, radius and surface gravity), the fol- 
lowing quantities are also obtained; mass-loss rate, current 
rotational velocity, current stellar mass, and surface abundances 
(He, B, C, N, O, Ne, Na, Si). For the simulation presented in 
this pap er, we use the LMC stellar evolution grid presented in 
iPaper J as input for our population synthesis models. 
In Sec. 13.11 we discuss the implementation method applied, 
and justify the initial distributions in Sec. 13.21 The selection 
effects applied to both the observed and modeled population are 
reviewed in Sec. 13.31 Note that some aspects of this selection 
effects are linked to the discussion in Sec. 13.21 Sec. 13.41 and l3. 51 
discuss stellar evolution models and rotational mixing. 

3.1. Population synthesis 

STARMAKER simulates a prescribed number of stars, with ini- 
tial mass (Mi), initial rotational velocity (^zams), and age (t) 
drawn from the corresponding distribution function (IMF, initial 
distribution of rotational velocities, and SFH). For a given initial 
triplet (Mi, uzams, and f), four evolutionary sequences from the 
grid of stellar evolution calculations are selected, and then used 
to interpolate the stellar parameters using a 3D linear interpola- 
tion algorithm. 

For a given Mi and ^zams, the interpolation routine first de- 
termines the two closest masses above and below Mi in the 
model grid. Then, for each of the two, it chooses the two mod- 
els with DzAMS closest to the required value (Fig. [1). We found 
that the algorithm works best when the four selected model se- 
quences form two pairs with equal masses. For example, to in- 
terpolate an I8M0 model, two 16 M© and two 19Mo models 
enter the interpolation. An interpolated track is show in Fig. [T] 
(bottom). Initial masses and initial velocities may be irregularly 
spaced in the model grid, but for every new mass introduced 
into the model grid, models of several rotational velocities are 
required to cover the dzams -parameter space considered in the 
simulation. 

As indicated by the arrows in Fig.lT] we first interpolate mod- 
els of equal mass to obtain models with the desired initial ve- 
locity, and then we interpolate in mass. We tested our result by 
changing the interpolation order, and found no significant dif- 
ferences. Adopting the nearest neighbors for interpolation , i.e 
the smallest distance -\/(AAf)2 -1- (A^zams)", would often result 
in selecting models of three or four different initial masses due 
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Fig. 1. Top: Illustration of the interpolation procedure of 
STARMAKER in the uzams - mi plane. Individual model se- 
quences are indicated by filled circles. Red lines indicate lines 
of constant mass. In this example, the properties of a star with 
initial mass Mi = I8M0 and initial velocity ^zams - 240km/s 
are determined through interpolation (the intersection of the blue 
dashed lines). STARMAKER chooses the circled models when 
interpolating. Models marked by triangles would be chosen if the 
nearest neighbor method was used. Black arrows indicate the se- 
quence of the interpolation, as described in the text. Dashed lines 
are intended as guide to the eye to distinguish models of similar 
velocity. Bottom: Interpolated track (red crosses), resulting from 
the example given above. The over plotted full blue line shows 
the model calculated with our stellar evolution code. The four 
black lines show the evolutionary tracks from which the model 
has been interpolated. 



to the domination of the velocity term in the distance equation 
(models marked by a triangle in Fig.[T). We found this could lead 
to unsatisfactory results. 

Stars with different initial masses and rotation velocities 
achieve different main-sequence lifetimes. To facilitate the in- 
terpolation in time, we use the fractional main-sequence life- 
time as the interpolation variable. For all models in the grid 
their MS-lifetimes have been tabulated so that, for any pair (Mi, 
fzAMs), the expected MS-lifetime can be interpolated. To en- 
sure a smooth behavior of isochrones in the HR diagram we 
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fit the MS-lifetime for each initial mass as function of lizams 
with a polynomial, and add a linear fit for masses where quasi- 
chemically homogeneous evolution is achieved for the highest 
rotational velocities (Fig.|2]i. The increase in the MS-lifetime for 
higher rotational velocity is due to rotational mixing of hydro- 
gen into the convective core. This effect saturates when the rota- 
tional mixing timescale becomes shorter than the main sequence 
timescale. At th is point the star becomes quasi-chem ically ho- 
mogeneous (e.g. lMaeded[T987l:lYoon & Lanpeil2005b . implying 
its interior is well mixed and the star can bum basically all the 
fuel it contains. 




2 3 4 

vrot/l 00 [km/s] 



5 6 



Fig. 2. The logarithm of the MS lifetime is plotted against the 
ZAMS rotational velocity for our 13 Mq (top) and 35 M© (bot- 
tom) sequences. In black we show the polynomial fit to the 
model data. The bottom panel shows also the linear part (dashed 
line) used to fit the models evolving quasi-chemically homoge- 
neously. 



Any observable can now be interpolated using the four 
model sequences used for interpolation (Fig. [TJ. As each model 
sequence is highly resolved in time (with several thousand time 
steps for the main sequence evolution), the interpolation error in 
time is negligible. 

3.2. Initial distribution functions 

In our population synthesis we draw the three initial parameters 
(Mi, uzAMS, and f), from the initial distribution functions, which 
are discussed below. We also adopt an inclination angle for the 
star, assuming random orientation in space. 

3.2.1. Star formation history 

We have used our population synthesis code to generate 
isochrones for non-rotating stars with ag es betwe en and 30 
Myrs from our stellar evolution grid (see iPaper lb . In order to 
evaluate the appropriate star formation history for our popula- 
tion synthesis model, we compare those isochrones with the dis- 
tribution of the stars of our observational sample in Fig. |3] 

The Nil region (blue symbols) contains a number of stars 
with ages younger than ~5 Myr These are primarily the massive 
O-type stars, associated with the recent star formation in, e.g.. 



LH 10. Both fields contain older populations with ages in the 
range of 5 to 25 Myr, and with no obvious preference for a par- 
ticular age (Fig.O. Given that we do not include the O-type stars 
in our analysis (see Sec. 13.3b . we therefore adopt continuous star 
formation as a good approximation of the observed lower-mass 
population in our population synthesis models. 
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Fig. 3. Isochrones generated from our non-rotating evolutionary 
models, compared to the location of stars of our sample, in the 
H-R diagram. Stars for which nitrogen abundances are available 
are shown with filled symbols, while we used empty symbols 
for those for which this is not the case. Black and blue symbols 
show data from NGC 2004 and Nil, respectively. 
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Fig. 4. Mass functions from the observed sample and our popula- 
tion synthesis simulation. The input initial mass function (IMF) 
used in the population synthesis is shown in green. The blue 
solid line represents the present day mass function (PDMF) of 
all stars in the simulation. The blue dashed line gives the PDMF 
of the stars which successfully pass our selected effect filters, in 
comparison the PDMF of all B-type stars of the observed LMC 
sample, including those for which no nitrogen measurement is 
available (red line) . The synthesized mass functions have been 
scaled by a factor of 1/500 for ease of comparison with the ob- 
served distribution. 
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3.2.2. Initial mass function 

As inp ut for our population synthesis we adopted a ISalpeterl 
(Il955h IMF. The masses are randomly picked from the IMF in 
the range between 7 and 60 M©, and the result is plotted in Fig.|4] 
(green, solid line). The present day mass function (PDMF) of our 
population synthesis simulation is also plotted (blue, solid line). 
It contains only main sequence stars and is calculated assuming 
continuous star formation. A comparison between the IMF and 
PDMF shows that mass loss has populated the 6 Mo-bin and that 
the (negative) slope of the PDMF has increased. The latter is an 
evolutionary effect: the most massive stars evolve more quickly 
and once they leave the MS they are removed from the popula- 
tion. At the high mass end the PDMF gets incomplete as mass 
loss shifts stars to lower mass bins and we do not have stars ini- 
tially more massive than 60 M© in our simulation. 

To test if our assumptions on the IMF and the SFH are rea- 
sonable, we construct the PDMF of the stars that would pass all 
observational selection effects in our population synthesis simu- 
lation (see Sec. l3.3l below). This distribution corresponds to the 
light-blue dashed line in Fig. |4] and is compared to the PDMF 
of the stars in our sample (red line). A plateau at lower masses 
develops, since many low mass stars are too faint to be included 
in the observational sample. At higher masses, the slope steep- 
ens further, since the hot O-type stars drop out of the sample 
(see Sec. 13.3b . This can be directly compared to the observed 
sample PDMF, shown as a red, dashed line in Fig.|4] Keeping in 
mind the low number statistics in the observed sample, it appears 
that below 30 M© the two curves have approximately the same 
shape, thus supporting the assumptions made for the simulation. 

The lowest mass stars in the observed sample are 9 M© while 
the lowest mass stars in the simulation are TMq. This is con- 
firmed in Fig. |6] where no stars are observed at the lowest 
masses. This could imply that the SFH is limited to ~25 Myrs, or 
— perhaps more likely — that the sample selection concentrated 
on the brighter stars in the field. This discrepancy between our 
model and the observational sample might introduce an error in 
our analysis. However, the error caused in the population num- 
bers of the different regions of the Hunter-diagram by including 
stars older than 25 Myrs in the simulation is on the order of 1 % . 
We discuss this in more detail in Sec. 14. II 



3.2.3. Initial velocity distribution 

iHunter et al.l (l2008bh fitted the observed rotational velocity dis- 
tribution of the LMC sample with a Gaussian function, exclud- 
ing radial velocity variables, supergiants and stars above 25 Mq. 
This mass threshold was used to avoid effects from mass loss 
induced spin-down on the velocity distribution (see Fig. [Sj. As 
such, all but the seven lowest mass O-type stars were not in- 
cluded in determining the velocity distribution. Be-stars below 
the mass threshold are i ncluded into the fit to the velocity distri- 
bution. As discussed b\l Hunter et al.' ('2009b|), the width of this 
function given by Hunter et al. (2008b) is too small and needs 
to be corrected by a factor (;r/4)'''*. Taking this into account, 
our ZAMS rotational velocities are selected randomly from the 
Gaussian function given in Eq.[T] (using /j. = lOOkm/s and cr - 
141 km/s ), which has been tru ncated at zero km/s and renormal- 
ized (see lHunter et al.ll2008bh . 



f(v) oc exp 



(1) 




150 200 250 
V sini [km/s] 



400 



Fig. 5. The projecte d rotational velocity distribution found by 
iHunter et al. I (l2008bl) for non-binary, non-supergiant stars of this 
sample with masses < 25 M© (black line). This is used as our 
initial rotational velocity distribution. The gray histogram repre- 
sents the rotational velocity distribution of our simulated pop- 
ulation, after application of all selection effects described in 
Sec. [331 



By using this function for our population synthesis, we as- 
sume that it is representative of the distribution function of the 



initial rotational velocities. This is justified when the stars rotate 
close to solid body rotation, which is ensured by including mag- 
netic fields in our models. Also, due to the low mass loss rates for 
the bulk of the stars in the sample, our models show little change 
of the rotational velocity over the main sequence evolution. To 
test this assumption, we have constructed a histogram of the 
present day projected rotational velocity distribution (PDVD) for 
the B-type stars from our population synthesis (gray bars in Fig. 
|5]l. Comparing this for the in itial projected rotatio nal velocity 
distribution (black line) from IHunter et al.l (l2008bb . our simu- 
lation predicts only slightly more slow rotators than the input 
function, but also slightly fewer fast rotators around 200 km/s. 
The overall shape at lower rotational velocities is preserved in 
our simulated data. Overall, the observed B-type star velocity 
distribution seems to be a good approximation to the initial ve- 
locity distribution of our population synthesis simulation. 



3.3. Selection effects 

The VLT-FLAMES survey was designed to contain a largely un- 
biased sample of O- and B-type stars. However, owing to both 
observational and analytical constraints, there are a number of 
selection effects present in the final Hunter et al. (2008b) sample 
(see|2]i, which we attempt to include in our simulation. Specific 
affected regions in the H-R diagram are highlighted in Fig.|6] 



3.3.1. Magnitude limit 

Our LMC sample is magnitude limited, to ensure sufficient 
signal-to-noise in each observed spectrum, with about V = 
15.53™ for the faintest star While the limit is not exactly the 
same in Nil and NGC2004, the limits in each are not signif- 
icantly different and hence it is reasonable to model these to- 
gether In our simulation, we therefore remove all simulated stars 
that are fainter than this threshold from our simulation, with this 
region marked in gray in Fig. |6] Magnitudes for the simulated 



5 



I. Brott et al.: Simulating a Population of LMC early B-type Stars as a Test of Rotational Mixing 




55 50 45 40 35 30 25 20 15 10 

Teff [kK] 

Fig. 6. H-R diagram of our population synthesis using the initial distributions discussed in Sec. 13.21 Regions subject to selection 
effects (Sec. 13.3b are highlighted as: gray - stars too faint for the observed magnitude cutoff; green - stars with Ten > 35 kK; red 
- stars with low surface gravities (log g < 3.2); yellow - remaining stars. Results from lHunter et al.l (l2009al) are over-plotted (blue: 
stars in Nil; black: NGC2004; with their shapes indicating their surface gravities). Radial velocity variables and Be-stars that are 
indicated by crosses and diamonds, respectively. Evolutionary tracks of lOM©, 20 M© and 30 of non-rotating models are shown 
for orientation. 



stars have been calculated using: 

V = M;„i+^ + Rv-E(B - V)-BC (2) 

M;;^! = -2.5 • log LModel + Mbol,0 (3) 

where Mj!j^j is the absolute bolometric magnitude of the sim- 
ulated star Observational parameters a re adopted s u ch tha t 
they are consistent with those used by iHunter et al.l (l2008bl) 
for NGC2004, i.e., distance modulus (//) = 18.56", E{B - 
y) = 0.09™, and a ratio of total-to-selective extinction Ry = 3.1. 
Bolometric corrections (BC) for Teff < 28172K have been 
adopted fro m Balona (11994 . and for Teff > 28172 from 
IVacca et all (Il996l) . 

3.3.2. Exclusion of 0-type stars 

The dominant nitrogen lines in O-type stars are those from N iii. 
However, this ion is notoriously difficult to model in non lo- 
cal thermodynamic equilibrium due to the need to include di- 
electronic recombination channels. To date, no large studies of 
nitrogen abundances in O-ty pe stars are ava i lable, althou gh He 
abundances were studied by iMokiem et al.l (l2006l l2007h . Due 



to the lack of nitrogen abundances for stars hotter than 35 kK 
we exclude such hot stars from the simulation (green group in 
Fig.|6l), corresponding to spectral types earlier than 09. 



3.3.3. Exclusion of low gravity objects 

The bul k of the stars have masses in the range of 10 to 20 M©. 
IPa per J argues that the surface gravity at the end of the main 
seq uence evolution for these objects is \ogg ^ 3.2 (see also 
Hu nter et al l l2008hh . While there are some observed stars with 
lower log g, their nature is unclear (IVink et al.ll20Tol) . We there- 
fore exclude post-MS objects, by rejecting all simulated stars 
with gravities below 3.2 dex (red regime in Fig. |6}. Note that 
this implies that we also reject the final part of the main sequence 
evolution of stars with masses above ~ 20 Mq, since the applied 
overshooting leads to an extension of the main sequence to lower 
gravities for these objects. We chose to do this because we favor 
a well defined selection criterion and because the most massive 
objects have only a low impact on our statistics. Most of the ob- 
served stars rejected are luminosity class I supergiants, but this 
group also contains three class II and two class III objects. 
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3.3.4. Exclusion of Be stars 

iHunter et al.l (l2009ah did not derive nitrogen abundances for 
most of the observed Be-type stars (the purple group in Fig. lA.lb . 
For consistency, we excluded a fraction of the most rapid ro- 
tators in our simulation. As the mechanism that produces the 
Be-phenomenon is still being discussed, it is difficult to define 
a clear criterion to identify simulated stars as Be stars. The Be 
mechanis m appears to b e linked with fast rotation and pulsa- 
tions (e.g. ICranmeijl2005h . Commonly quoted required fractions 
of critical ro tation for the occurrenc e of the Be phenomenon 
are 70-80% ("Porter & Riviniusll2003h . though values as low as 
> 40% (Cranmer 2005) are also d erived. iMartayan et al. (2008) 
finds Q/Qc = 85% in the LMC, and lTownsend et al.l ( f2004 argue 
that Be-stars might actually rotate as fast as 95% of the breakup 
velocity, since the line broadening due to rotation becomes in- 
sensitive to gravity darkening at high rotational velocities. 

The fraction of Be/(B+Be) in the observed sample is 13.5%. 
We can define a theoretical Be fraction within the limits of the 
yellow region of Fig. |6]by dividing the number of stars rotat- 
ing faster than a specified fraction of break-up rotation by the 
number of all stars in that region. To match this fraction with 
the observed one in our sample would require a velocity limit of 
55% of the break-up velocity. This number seems too low. As 
clearly more effects than fast rotation are involved in producing 
Be stars, we decided to not try to reproduce their number in the 
population synthesis. In an attempt to simulate the effect of the 
exclusion of Be-stars from the sample, we chose a more conser- 
vative limit and excluded stars with a rotational velocity higher 
than 90% of their breakup velocity in our simulations. This cor- 
responds to 1 .6% of all stars in the yellow region of Fig. |6] 

In summary, we are unable to reproduce the fraction of ob- 
served Be-stars using our cut-off of 90% of breakup velocity. 
However as we only exclude a small minority of stars from the 
simulation, we do not expect that this will be a significant source 
of error. 



3.4. Stellar evolution models 



We calculated a grid of 284 stellar evolution models past 
core hydrogen exhaustion using a one-dimensional hydrody- 
namic stellar evolution code, that includes the ph ysics of dif- 
fere ntial rotation and m ag netic fields as describe d in lPaper j : see 
also lYoon et"al1 (l2006h and lde Mink et al.l(l2009l) . The grid spans 
initial masses of 5-60 Mq and rotational velocities (wzams) of 0- 
600km/s, as shown in Fig.|7] The grid is somewhat irregular in 
rotational velocity since we initialize our models as rigid rota- 
tors in a chemically homogeneous state of thermal equilibrium. 
The subsequent adjustments to achieve CN-equilibrium in the 
core alter the stellar structure on a thermal time scale and lead 
to a somewhat different surface rotation rate, which thereafter 
changes only on the nuclear time scale, and which we designate 
as uzAMS- 

The effect of the centrifugal force on th e stellar structure 
is considered following Endal & SofiiJ (Il976l) . We consider the 
transport of angular momentum and chemical elements due to 
various rotationally induced hydrodynamic instabilities, which 
include the Eddington-Sweet circulation, the dynamical and sec- 
ular s hear inst ability, and the Goldreich-Schubert-Fricke insta- 
bility (iHeger et al. 20Q(^. Th e transport o f angular momentum 
by the Spruit- Tayler dynamo (ISpruitll2002l) is also implemented. 
Angular momentum transport and chemical mixing is approxi- 
mated as a diffusive process. 

Convection is treated by the mixing length theory assum- 
ing a m ixing-length parameter of umlt =15. We follow lEanged 
(Il99lh for semi-convection, for which rather fast mixing is as- 
sumed, with a semi-convection efficiency parameter of ffsEM^l- 
We consider convective overshooting using an overshooting pa- 
rameter of 0.335 Hp. This value results from our new calibration 
using the observed u sin / drop that is foun d in our data when we 
plot V sin ; against the surface gravity (see iPaper II) for more de- 
tails). The rotational mixing efficiency has been calibrated such 
that the trends in nitrogen of our sample are reproduced well. 
This leads to adopting /c =0.0228, which is the ratio of the 
turbulent viscosity to the diffusion coefficient (see iHeger et aH 
120001) . 

We use the recipe of IVink et all (l2000i l200l for mass loss 
via stellar winds from O- and B-type stars. To account for the 
effects of h elium enrichment o n th e stellar wind, we follow the 
approach of lYoon et al.l (l2006l) (see iPaper II for details). 
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Fig. 7. Parameter space covered by the LMC evolution grid. 
Each model sequence is represented by a dot marking its initial 
mass and rotational velocity. 



3.4.1. Initial chemical composition 

We choose an LMC mixture for the initial chemical composition 
of our models ( Table[l]l. CNO abunda nces are taken from results 
for H II regions dKurt & Dufou]J'1998). Mg and Si measurements 
fo r B-type stars from H unteret^al. (2008b) and B-supergiants 
of iTrundle et"aD ilOOlh are, on average , 0.37 dex lower than 



the solar values of lAsplund et al. ( 20051). Therefo re, we have 



scaled all other elements from Asplund et al.l (l2005 j ) bv -0.4 dex 
and adopted the average values of Mg and Si from lHunter et al.l 
(2007|) andlTYundle et al. (2007) for our stellar models. The re- 
sulting metallicity of this mixtur e is Z = 0.0047. We have 
adopted the OPAL opacity tables (llglesias & Rogersl[T996l) . for 
which we use the Fe abundance to interpolate between tables of 
different metallicities. 

We have adopted the helium mass fraction using the pri- 
mordi al helium abundance (Fp= 0.2477) from iPeimbert et all 
(I2007h . Assuming the helium abundance to be a linear func- 
tion of the metallicity Z and fixing it to 0.28 at solar metallicity 
(iGrevesse et al.llI996l) results in a value of 0.2562 for our LMC 
models. 
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Table 1. LMC mixture used in our evolutionary models. 



c 


7.75 


Mg 


7.05 


X 


0.7391 


N 


6.90 


Si 


7.20 


Y 


0.2562 





8.35 


Fe 


7.05 


Z 


0.0047 



Notes. Tlie baseline abundances of C, N, O, Mg, Si and Fe are given 
in the first two columns, i.e. log(C/H) + 12. The total hydrogen (X), 
helium(Y) and metals (Z) mass fractions are given in the third column. 



More details on the calibration and choice of the chemical 
composition for the full evolutionary model grid, which also 
include s models for Galactic and SMC stars, are presented in 
iPaperJ . 



3.5. Effects of rotational mixing 




5 10 15 20 25 

Age [Myr] 

Fig. 8. Nitrogen surface abundance as a function of time for 
models with different rotational velocities for 10 M© (green), 
13 Mo (blue), 20 Mq (black) and 40 Mq (red) are shown. The 
model calculations have been ended shortly after the end of core 
hydrogen burning. The ZAMS velocities are indicated next to 
each track. 

In our models, the strong transport of angular momentum 
due to the Spruit-Tayler dynamo keeps the star rotating al- 
most rigidly on the main sequence. As a consequence, the 
she ar instability - which is iinporta nt for non-magnetic mod- 
els (iHeger et al.ll200(it lMaedei]l2000l) - is less efficient than the 
Eddington-Sweet circulation in our models. 
The surface abundances of stars can change due to rotational 
mixing on the main sequence but not all elements are equally 
well suited for tracing such mixing. While helium is produced 
in large amounts in the core, the models for moderate rotational 
velocity show only small enhancements of helium at the surface. 
This is because the chemical stratification built up by hydrogen 
burning between the helium enriched core and the hydrogen-rich 
envelope effectively hinders chemical mixing across the bound- 
ary. Only the fastest rotating models (dzams > 350km/s) be- 
come significantly helium eniiched at the surface, and some 
of them even follow the quasi-chemically homogeneous evolu- 
tion if rotational m ixing operates faster than hydrogen burning 
(lYoon et al.ll2006h . Additionally, helium is already abundant at 



the stellar surface initially, and the expected enhancement rela- 
tive to the initial helium abundance is typically small. 

Nitrogen, however, can be considerably enhanced at the stel- 
lar surface on the MS, even in relatively slow rotators. At the 
same time, carbon is depleted, but with a smaller relative change 
than nitrogen, and oxygen remains almost constant. This can 
be understood in the following way. At the onset of the CNO- 
cycle, there is a shoit phase of the CN-cycle in which most of 
the carbon is transformed to nitrogen. Duiing this phase, hardly 
any helium is produced, and the mean molecular weight gradi- 
ent at the outer edge of the core remains small. Therefore, ni- 
trogen produced by the CN-cycle can be efficiently transported 
into the layers above the convective core. Once the mean molec- 
ular weight gradient increases due to helium production, such 
mixing becomes much slower As mixing throughout the radia- 
tive envelope — i.e. above the major mean molecular weight 
barrier produced by hydrogen burning — continues during the 
main sequence evolution, the surface becomes gradually en- 
riched with nitrogen. As the nitrogen abundance at the surface 
of early B type stars can be measured from optical spectroscopy, 
nitrogen is one of the best tracers for rotationally-induced chem- 
ical mixing. 

In Fig. [8] we show the surface nitrogen abundance as a func- 
tion of time for a selection of models from our grid. As the 
nitrogen eniichment during the main sequence evolution is a 
monotonic function of time, nitrogen might be used as an age 
indicator for stars with known initial stellar parameters. In our 
models, the nitrogen abundance achieved at the end of core hy- 
drogen burning is higher for higher initial rotation velocities, as 
the Eddington-Sweet c irculation operates at a higher flow veloc- 
ity for faster rotation (lYoon et al.ll2006h . At the same time, at 
fixed initial rotation velocity, the maximum nitrogen abundance 
increases only slightly with initial mass in the range of 10 to 
20 Mq, but it becomes more sensitive at masses above 30 Mq. 
When hydrogen burning stops, the core contracts and the enve- 
lope expands and eventually becomes convective, leading to a 
further dredge-up of nitrogen and a sharp increase of the nitro- 
gen surface abundance as shown in Fig. [8] 

The abundances of other elements also show variations due 
to rotational mixing. Sodium is enhanced due to proton captures 
in the NaNe-cycle, but it is difficult to determine observationally 
as its lines are weak in this mass range and they are contaminated 
by interstellar absorption. Boron is strongly depleted, because it 
is destroyed as it is mixed fr om the r elatively cool outer layers 
into hotter inner layers (see iPaper Ih . but its spectral lines are 
weak due to the low initial abundances. Carbon depletion can 
also be observed in the optical, but it ca nnot be measured a s 
accurately as nitrogen due to NLTE effects ( iHunter et al.l2008"bl) . 
Fluorin e and alu minum are also affected. A full account of this is 
given in lPaper j : in what follows we concentrate on the nitrogen 
abundances. 



4. Results 

To investigate the role of rotational mixing in massive stars, 
Fig.|9]shows projected equatorial rotational velocities versus ni- 
trogen abundance (hereafter the Hunterplot). The figure contains 
both observations and theoretical predictions. The observational 
data includes all stars from our LMC sample for which nitrogen 
abundances are available (see Sec. |2]i. Radial velocity variables 
are marked by triangles, while stars with no detected variations 
are represented with circles. The color coding of the symbol in- 
dicates the surface gravity of each object as shown. 
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3.2 < logg < 3.7dex • logg>4.1dex • 

3.7 < logg < 4.1dex • 

Fig. 9. Hunter-plot showing projected rotational velocity against nitrogen enhancement. Our population synthesis is sh own as a 
density plot in the background. The color coding corresponds to the number of stars per pixel. Overplotted are data fromHu nter et al.l 
(l2009ah . with surface gravities as indicated by the color (see figure key). Single stars are plotted as circles, radial velocity variables 
as triangles. Evolutionary tracks of 13 Mq, corresponding to the average mass of the sample stars, are shown with their surface 
gravity coded by the same colors as the observations. The velocities of the tracks have been multiplied by 7r/4, to account for the 
average projection effect. The cross in the lower right corner shows the typical error on the observations. 



The simulation considered 1.5 million stars with initial 
masses between 7 and 60 M©. Only those that remained in the 
main sequence can be compared to the observational sample. 
These are plotted in Fig. |9] About 45.7% of initial sample have 
left the main sequence, whilst another 46.3% were excluded due 
to selection criteria based on magnitude, temperature (spectral 
type) and surface gravity (see Sec. 13.3) . The remaining sample 
contains 120728 stars. 

The most important rejection criterion is the magnitude cut- 
off, but the age of the simulated star also has a considerable ef- 
fect. This is easily understood considering that the maximum age 
of 35 Myr in our SFH corresponds to the MS-lifetime of an 8 
model. As the MS-lifetime decreases rapidly with increasing ini- 
tial mass (see Fig.©, the chance that a simulated star is assigned 
an age greater than its MS-lifetime becomes significant. 

The magnitude cutoff eliminates a large number of simulated 
stars because the IMF favors lower masses. For instance, models 
with less than 20 M© are too faint to be observed for most of 
their MS-lifetime. A I2M0 star spends 72% of its MS-lifetime 
below the magnitude limit, while a 15 M© star is not observable 
for 35% of its MS-Hfetime. 

The simulated data have been grouped into 5 km/s x 0.04 dex 
bins to generate a density plot. The color coding represents the 
number of stars in the bin (see the color bar on the right in Fig.|9]l. 
We also show evolutionary tracks of 13 M© for comparison. For 
these tracks the surface gravity is indicated with the same color 
coding as for the observational data. 

Fig. [To] shows the same population synthesis simulation as 
Fig. |9] but we have added a random Gaussian error with cr = 



0.2 dex to the predicted nitrogen abundances. This error is char- 
acteristic of the spread in the nitrogen abundance determinations 
of the unenriched stars. Note that it is smaller th an the mean 
nitrogen erro r of the entire sample (0.3-0.4 dex; iHunter et al.l 
'2OO7I l2008bh . which is represented by the cross in the lower 
right corner of Fig. |9]and [TO] We did not include an error for 
the predicted u sin / values as the observational values should be 
reasonably secure. The color coding of the observational data in 
Fig.[TO]represents their estimated evolutionary mass. 
We have divided the diagram into five regimes: 

Box l(ab)isa region of the diagram at high u sin / which is pre- 
dicted to be almost empty. The full line defining box l(ab) 
gives the locus where the number of simulated stars per pixel 
is ~ 40, while this number is ~ 15 for the dashed Une defin- 
ing Box lb. 

Box 2 is a region of the diagram where very few stars are pre- 
dicted, but where a group of enhanced, slowly rotating stars 
are observed . This g roup of stars was already described by 
iHunter et al.l (l2009ah . and will be discussed further below. 

Box 3 encloses the stars in our simulation which have a clear 
signal of nitrogen enhancement due to rotational mixing. 

Box 4 is a region where neither stars are observed, nor a sig- 
nificant number of them is predicted. 

Box 5 is the region of the diagram where the un-enriched stars 
are located. The initial nitrogen abundance is about 6.9 with 
a typical spread of about 0.2 dex. We chose the upper bound- 
ary of 7.2 such that the probability of observed stars found 
above this line to be truly enriched is large. 
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Fig. 10. As for Fig. |9]but with a random error (cr = 0.2 dex) added to the nitrogen abundances in the simulated data. The observa- 
tional data have been color coded according to their evolutionary masses (see figure key). For the boxes indicated, in Table |2] we 
give the ratios of the number of observed to simulated stars 



Table 2. Fraction (in per cent) and absolute number of observed 
and simulated stars found in each box defined in Fig.fTOl 







observed 




modeled 




all 




excl. UL 








% 


# 


% 


# 


% 


#/103 


Box lab 


14.0 


15 


7.4 


6 


4 6+"" 


5.6 


Box la 


5.6 


9 








1 6+"^ 


2.0 


Box lb 


8.4 


6 


7.4 


6 


3 0+"^ 


3.6 


Box 2 


15.0 


16 


19.8 


16 


1 9+0.0 
^•^-0.1 


1.4 


Box 3 


17.8 


19 


13.6 


11 


35 7+-'-'' 

■^-'•'-4.6 


43.1 


Box 4 












3.1 


Box 5 


53.3 


57 


59.3 


48 


55 T'-^ 

-'-'•'+6.2 


67.5 


Total 


107 


81 


120.7 



Notes. For the observed stars, in Column 2 we take all stars into ac- 
count, in Column 3 we only consider those without upper limits. The 
third column gives the expected number of stars from our population 
synthesis model. Errors on the number of simulated (Column 6) stars 
are derived by comparing to simulations with a velocity distribution 
that was broader or narrower by ±25 km/s. 



In Table |2] we present the absolute and relative number of 
stars populating each box in both the observations and simula- 
tions. For the observations we present statistics including and 
excluding objects for which only an upper limit to the nitrogen 
abundance was established. A downward revision of the upper 
limits might reduce the number of observed stars in Box 3, while 
that in Boxes 5 and 1 might increase. For all stars in Box la we 
only have upper limits, but their presence in this box is not af- 
fected by this. Unless otherwise stated, the numbers quoted in 



the discussion below include the stars with upper limits to their 
nitrogen abundances. 

We find 53.3% of the total observed sample to be in Box 5, 
containing the unenriched stars which do not rotate rapidly. This 
is in good agreement with the theoretical prediction of 55.9%. 
Due to the large number of stars in this region, this agreement 
is not affected by stars with upper limits on their nitrogen abun- 
dance. Most of the stars in this box are expected to have a rel- 
atively low mass, and in an advanced phase of core hydrogen 
burning such that they have become sufficiently bright in the op- 
tical to be above the magnitude cutoff. Indeed the observed stars 
with M < 10 Mq accumulate in Box 5 and appear to have mainly 
surface gravities of less than 3.7. 

We note that a more efficient mixing would reduce the pre- 
dicted number of stars in Box 5. However, the mixing efficiency 
is calibrated to best represent the correlation between nitrogen 
abundance and v sin / (see above). The fact that after such a cali- 
bration the number ofnon enriched stars predicted in Box 5 is in 
very good agreement with the observations supports the theory 
of rotational mixing. 

Box 3, containing the stars of our simulation which show 
rotationally induced nitrogen enhancement, is predicted to con- 
tain 35.7% of our sample. However only 17.8% of the observed 
sample is found in this region. Hence we predict approximately 
twice as many stars to be enriched as are actually observed. 
Subtracting the stars for which we have only upper limits leaves 
only 13.6% of the sample stars in Box 3, increasing the discrep- 
ancy between theory and observations to almost a factor of three. 
On the other hand. Fig. [TT] shows that more observed stars in 
Box 3 are found in the \2...2QMq bin compared to the lower 
mass bin (M < 12 M©), while our models predict similar (abso- 
lute) numbers for both mass bins. This may indicate that subdi- 
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viding our sample into too small subsamples results in a loss of 
statistical significance. 

While this may show that analyzing subsamples may lead to 
effects of small number statistics, it also shows that the overall 
discrepancy in the predicted and observed number of stars in 
Box 3 may not be very significant. 

Box 2 contains a group of apparently slowly rotating but 
nitrogen enhanced stars. Theoretical models of rotating single 
stars predict such enhancements only for fast intrinsic rotation 
(fZAMS > 150km/s). Note that all the observed targets in Box 
2 cannot be rapid rotators seen nearly pole on. Indeed our sim- 
ulations (which assume random inclination angles) predict that 
only 1.2% of the sample (corresponding to one star) would be 
found in this region. Observationally 16 stars comprising 15% 
of the sample lie in Box 2. Therefore the projection effect can- 
not populate this part of the diagram. This statement is enforced 
by the fact that no stars are observed in Box 4, which — if the 
projection effect were to explain the large number of stars in 
Box 2 — should contain more stars than found in Box 2. We 
conclude that the vast majority of stars in Box 2 are intrinsically 
slow rotators and their high nitrogen abundance is not explained 
by our models. The nitrogen enhanced stars in Box 2 appear to 
have higher masses and to be less evolved. There are no stars 
with M < IOM0 in this group, even though this is one of the 
most populated mass bins in the observed sample. With the ex- 
ception of five stars (of which one is a binary) surface gravities 
in this group are larger than 3.7 dex. A comparison with the evo- 
lutionary tracks in Fig.|9]indicates that most stars of this box are 
in an intermediate or evolved stage of their MS-evolution. The 
large number of stars found in Box 2 is not predicted by the the- 
ory of rotational mixing in single stars with an overpopulation of 
observed targets by more than a factor often with that predicted. 

Box l(aH-b) contains a group of 15 fast rotating but unen- 
riched stars (14% of the observed sample). We predict a factor 
of 3 less stars (4.6%) in this box, with most of them expected 
to be in Box lb, i.e. very close to Boxes 3 and 5. In particular 
the area where the six observed stars in Box la are found is pre- 
dicted to be empty. This prediction is largely due to the adopted 
observational selection effects, but to a small extent also due to 
the fast nitrogen enrichment of rapidly rotating stars (Fig. [S). 
Most of the simulated targets populating the diagram started 
their evolution below the magnitude cutoff. Those which were 
rapidly rotating have sufficient time (before passing the magni- 
tude limit) to enhance their surface nitrogen abundance so sig- 
nificantly that once they appear in the diagram they show up in 
Box 3. In Box l(a or a+b), there does not seem to be a prefer- 
ence for high surfac e gravities and low masses as suggested by 
iMaeder et al.l (l2009l) . The six objects in Box la have masses in 
the range of 9 ... 15 M©, which is a typical spread for the sam- 
ple. Their average mass is 12 which is similar to the average 
for all the stars for which nitrogen abundances have been deter- 
mined of 13.6 Mo. Additionally the stars cover a range of ages. 
As the stars in Box lb have upper limits to their nitrogen abun- 
dances, some of them could actually lie in Box la, increasing the 
discrepancy with prediction. In summary the existence of a sig- 
nificant number of targets in Box la cannot be reconciled with 
the prediction of rotational mixing for single stars that this area 
should contain effectively no objects. 

Fig. [TT] indicates that whilst the population of Box lb is 
strongly mass dependent, that dependency is less pronounced 
for Box la. The top panel shows simulation and data with 
M < 12 Mq. The models of 12 Mq spend about 70% of their 
MS at luminosities below the magnitude limit. When they fi- 
nally become bright enough, they are already so enriched that 



they show up in Box 3. Only for models of intermediate mass, 
as shown in the middle panel (I2-2OM0), does the interplay be- 
tween the selection effects and the mixing-timescale allow stars 
to populate Box 1. Using an even finer mass binning, we find 
that Box la is actually only populated by stars in the mass range 
12-15 M0, while Box lb contains models in the mass range 12- 
20 Mq. The bottom panel shows models and observed stars more 
massive than 20 M©. In both cases Box 1 remains empty. Models 
of 20 M0 or more are above the magnitude limit all time, but they 
have effective temperatures above 35 000 K during their early 
main sequence evolution, such that the fast rotators, once they 
have become cool enough to appear in our simulation data, have 
akeady achieved a significant nitrogen enrichment. 

4. 1 . Model uncertainties 

Before we focus on the possible evolutionary implication of our 
findings, we first discuss the significance of the results summa- 
rized in Table |2] which are obtained from our population synthe- 
sis. The discrepancy between observation and prediction is very 
large for both Boxes 1 and 2. In contrast, it is smaller for Box 3, 
being only a factor of approximately two. Hence in discussing 
our uncertainties, we will focus on how they might affect our 
simulation of the fractional population of this Box. 

These numbers are subject to various uncertainties: i) the 
mixing efficiency in the stellar models, ;/) the initial velocity dis- 
tribution, mass function and star formation history. Hi) the obser- 
vational biases, and iv) exclusion of unsuitable stars. 

4.1 .1 . Mixing efficiency 

The cal ibration of the mixing efficiency is described in detail in 
iPaper ll The models in our stellar evolution grid are calculated 
with a mixing efficiency parameter of f- - 2.28 x 10"^ (see 
Sec. 13.4b . The impact of a larger mixing parameter f- is that a 
larger maximum nitrogen enhancement is achieved for a given 
initial rotation rate. Therefore, the general trend of nitrogen en- 
hancement versus v sin ; would be shifted to lower rotation veloc- 
ities and larger nitrogen abundances. In Fig. |9] the evolutionary 
tracks of our 13Mo models show an unaffected nitrogen abun- 
dance for u sin i < 100 km/s within the observational error, and a 
peak value of 7.9 is obtained at v sin ; ^ 300 km/s. This behavior 
does reflect the general trend of the observed stars in Boxes 5 
and 3 quite well: the maximum observed nitrogen abundance, 
the slope of the nitrogen versus u sin /-relation, and the foot-point 
at Dsin/ ^ 100 km/s are all matched well. We thus believe that 
the calibration of the mixing efficiency of our models does not 
contribute significantly to the observed discrepancy. 

4.1 .2. Input distribution functions 

The adopted distribution of initial rotational velocities has been 
discussed in Sec. 13.2.31 In order to assess the sensitivity of the 
predictions on the distribution of initial rotational velocities, we 
investigated the effect of a broader (cr - 166 km/s) and narrower 
(cr = 116 km/s) width of the underlying Gaussian . A variation 
of +25 km/s is consistent with the errors given in iDufton et al.l 
(l2006h on the fit to the projected rotational velocity distribution. 
The relative changes are given as errors on the theoretical per- 
centages in Table |2] While we find that the numbers predicted 
for most boxes do not change significantly, the narrower veloc- 
ity distribution leading to fewer rapid rotators, shifts about 6% of 
all modeled stars from Box 3 to Box 5. This reduces the discrep- 
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ancy for Box 3 (35.7% predicted vs. 17.8% observed) to below a 
factor of two (29.8% predicted vs. 17.8% observed), but it can- 
not eliminate it. 

The initial mass function is assumed to be Salpeter [a - 
-2.35; see Sec. 13.2.2b . Given the apparent universal nature of 
this distribution in the mass range from 10-60 Mo (IChabrieii 
l2003HKrouDa> 200lV and given the good agreement between the 
observed and simulated mass functions (see Fig. |4|i we beUeve 
that this is unlikely to be a significant source of error 

As discussed in Sec. 13.2.11 we have assumed continuous 
star formation in our population synthesis model. Clearly, the 
star formation history may influence the distribution of stars in 
Fig. [To] For instance, a young co-coeval population will cause 
the fraction of predicted stars in Box 5 to increase, at the ex- 
pense of the number of stars in Box 3. However, as Fig. |3] ap- 
pears largely consistent with our assumption of continuous star 
formation, we do not expect a significant error due to this as- 
sumption. We have investigated the case where star formation 
stopped after 25 Myr (as suggested in Sec 13.2) instead of contin- 
uing for 35 Myr The stars that are older than 25 Myr are mainly 
to be found in Boxes 3 and 5. If we exclude these stars, the num- 
ber of stars in Boxes 3 and 5 change by about 1%, and by much 
less in the other boxes. 

4.1.3. Observational biases 

While we accounted for the observational biases in our popu- 
lation synthesis model, some error may be associated with this. 
Most relevant in this respect may be the magnitude cut-ofF ap- 
pUed in the VLT-FLAMES Survey. While Fig. |6] shows that in- 
deed no observed stars fainter than our adopted magnitude cut- 
ofF are found in the sample, there may be a slight deficit of ob- 
served stars near the magnitude cut-ofF, especially at masses be- 
low -10 Mq. As these stars would be nearly at the end of their 
main sequence evolution, the ratio of enriched to unenriched 
stars amongst them would be greater than average. However, 
limiting the population synthesis to stars above 9 M© has a very 
similar effect as limiting the star formation history to 25 Myrs 
(discussed above), since the main-sequence lifetime of a non- 
rotating 9 Mq model is about 25 Myrs. It changes the number 
of stars in Box 3 and 5 only by less than one percent. The 
main-sequence lifetime of a non-rotating 9 M© model is about 
25 Myrs. 

Binary stars are slightly brighter than single stars of a mass 
or brightness equal to the binary primary. Therefore, some bina- 
ries with primaries which would not pass the magnitude limit by 
themselves might, through the additional light of their compan- 
ion. A test population simulation where we assumed all stars to 
consist of binary stars with identical companions convinced us 
that this effect does not significantly affect any results described 
below. 

For very rapid rotators, the observed brightness and gravity 
may also depend on the inclination angle due to the latitudinal 
dependence of gravity darkening. However, the observed and the 
modeled number of very r apid rotators is so sm all that this effect 
is negligible here (see also lHunter et al.ll2009ah . 
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Fig. 11. Mass dependence of the Hunter-plot. We have separated 
Fig.[TO]into three mass bins. Top to bottom: Population synthesis 
results and observational data for the mass bins M < 12 M©, 
12 < M < 20 Mq, M > 20Mq are shown, respectively. 



4 1 4 Removal of unsuitable stars ^^^y companion (e.g. I Sana et afll2009h which will modify the 

evolution during core hydrogen burning. Other m assive stars 
In our population synthesis simulation, we have only included are clearly affected by internal magnetic fields (e.g. IWade et aH 
single stars computed with a single set of physics assumptions. 2006). It remains unclear wether the magnetic fields used in our 
In reality however, different physical situations may occur For evolutionary models for angular momentum transport describe 
example, many massive stars are known to have a close bi- those detected in massive stars. Our population synthesis model 
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is therefore only representative of at most a certain fraction of 
the observed stars. 

In order to assess the potential consequences for our conclu- 
sions concerning the number of stars in Box 3, let us first assume 
that the stars in Boxes 1 and 2, which cannot be explained well 
by our model anyway, need to be discarded as "unsuitable" for a 
comparison with our model. This would reduce the total number 
of stars to be compared with from 107 to 76, and then render the 
fraction of observed stars in Box 3 to 25% (rather than 17.8%), 
which is closer to the predicted fraction of 35.7%. 

However, the fraction of observed stars which is "unsuitable" 
for a comparison with our model may be greater; let us assume 
it is as high as 50%. As the observed stars in Box 3 agree well 
with the physics of rotational mixing, let us further assume that 
none of the 50 or so observed stars to be disregarded would lie 
in Box 3. Again, we would assume that the observed stars in 
Boxes 1 and 2 should be disregarded. However, about 20 stars 
need to be disregarded also from Box 5. For example, this could 
be fast rotators which cannot mix because of strong internal 
magnetic fields. With such assumptions, we had 37 "suitable" 
stars left in Box 5, and 19 in Box 3. The latter number corre- 
sponds to ~33% of the reduced sample, which then is close to 
the predicted fraction of 35.7%. While the assumption that we 
can describe only half of the stars in the observed sample may 
be extreme, this case cannot be ruled out, and it serves as an ex- 
ample of the difficulty in defining a robust observational sample. 

4.1.5. Stochastic effects 

Here, we assess the spread in the number of stars obtained for 
each box in Fig. [10] if only 107 stars are drawn randomly from 
the simulated population. This is equivalent to the following 
experiment. A pot contains N=120728 balls, equivalent to the 
number of simulated stars in Fig. [10] The balls in the pot have 
the colors Ci . . . C5. The number of balls in the pot with color c, 
is equal to the number of simulated stars in Box i, for ; = 1 ... 5. 
We draw n=107 balls from this pot, without taking the order into 
account in which they are drawn. Because of the large number 
of balls in the pot, the probability p, to pick a certain color can 
be approximated as constant, i.e. as the number of balls of color 
c, divided by the total number of balls N. This approximation is 
valid as long the number of balls of color c, is much larger than 
the number of balls drawn from the pot. 

When the experiment i s repeated man y times, the result is 
multinominally distributed dPapoulisll 1 984 p. 75). The expecta- 
tion values /Y, for the number balls drawn from each color and its 
standard deviation cr, can be calculated by: 

/^i = npi (4) 
oj = npi{l-pi) (5) 

The expectation values given in Table [3] agree with the percent- 
ages found in Table [2] The standard deviation indicates a rather 
small spread, which cannot cause the discrepancy between ob- 
servation and simulation discussed above. This is true in partic- 
ular for the ratio of the number of stars from Box 3 and Box 5. 

4.1.6. Overall assessment 

We have shown above that the discrepancy of a factor of two in 
the predicted and observed fraction of stars in Box 3 is not suf- 
ficient to rule out rotational mixing for this subsample. Our the- 
oretical mixing calibration, the adopted star formation history, 
and initial distribution functions for mass and rotational velocity 



Table 3. Expectation value of the numbers of stars and standard 
deviation for each box 





fi ± o-m 


yu ± tr[%] 


Box 1 


4.9 + 2.2 


4.6 ± 2.0 


Box 2 


1.3 + 1.1 


1.2 + 1.0 


Box 3 


38.2 + 5.0 


35.7+4.6 


Box 4 


2.8 + 1.6 


2.6 ± 1.5 


Box 5 


59.8 + 5.1 


55.9 ±4.8 



Notes. Expectation values and standard deviations as calculated from 
equations [4] and |5] are given in column one. In the second column these 
values are expressed in percent of the sample. 

apparently induce little uncertainty. However, the adopted mag- 
nitude cut-off and in particular the possibility of applying dif- 
ferent physical situations to main sequence stars than assumed 
in our model allows for a modification of the fractions of stars 
in Boxes 3 and 5 which may bring theory and observations into 
agreement. 

On the other hand, the possible agreement of the predicted 
distribution of stars in Boxes 3 and 5 is not sufficient to verify 
the theory of rotational mixing. This will be discussed further 
when predictions of massive close binary evolution are included 
in our considerations, in the next section. 

5. Discussion 

From the previous section, we conclude that overall, the agree- 
ment between the population synthesis model and the observed 
sample is rather poor While the model predicts stars to be con- 
fined to Boxes 5 and 3, the observed stars are much more wide 
spread. Quantitatively, only in Box 5 is the number of observed 
stars in good agreement with the predictions. While the stars in 
Box 3 at least qualitatively agree with the population synthesis 
model, the stars in Box 1 and Box 2 are completely unexpected 
and cannot be explained by the model. Before any further con- 
clusion, it is thus worthwhile to consider physical and evolu- 
tionary processes which could populate Boxes 1 and 2, as those 
might have effects on Boxes 3 and 4. 

5.1. Close binary evolution 

Though our model assumption is that the diagram is popu- 
lated only by single stars, it is well known that a significant 
fraction of stars is part of a binary system (e.g., iMason et al.l 
2009), with a bi nary fraction of 0.5 or more found in y oung 
star clusters (e.g.. lSana et al.l2008l [20091: iBosch et al.l2009h . and 
with a period distribution such that many of them will interact 
dSana & Le Bouquinl2010l) . We emphasize that strong binary in- 
teraction (i.e., mass transfer or merging) will often lead to a sys- 
tem that is either a single star or which shows insignificant ra- 
dial velocity variations. Most observed radial velocity variables 
in the VLT- FLAMES Survey will likely be pre-interaction sys- 
tems and thus well suited for comparison with single star evo- 
lution. The apparent binary fraction in our various observational 
subsamples thus give no direct clue to the importance of strong 
binary interaction for the subsample. 

In close binaries, the primary star may fill its Roche lobe 
during the main sequence evolution. This results in the removal 
of the stellar envelope of the primary star, which is transformed 
into a (mostly unobservable) helium star If the entire envelope 
is added to the companion, the surface nitrogen abundance of 
this star will increase by a factor of 3 to 5, even when rota- 
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tional mixing is not taken into account (cf.. Fig. 9 in lLanger et al.l 
I2008h as in the accretion process the envelope material is well 
mixed. Subsequent rotational mixing of the mass gainer may en- 
hance the surface nitrogen a bundance further on the thermonu- 
clear timescale (cf.. Fig. 10 in lLanger et al.l2008b . If only a small 
fraction of the material is accreted, the result is no or only a mild 
N enrichment. During the accretion process not only mass but 
also angular momentum is transferred to the companion, and it 
has been found that the accretion of only about 10 % of the stel- 
lar mass is able to spin-up the star significantly (|Packetlll98lt 
iLanger et al.ll2008l) . 

Therefore, conservative and moderately non-conservative 
mass transfer would tend to populate the upper right corner 
in Fig. [T Ol (Box 3), whi l e high ly non-conservative mass trans- 
fer (e.g.. IPetrovic et all I20p5[) is likely to m ove stars to the 
lower right corner (Box 1) ( Langer et al.ll2008l ). Binary interac- 
tion leading to a merger is also expected to result in fast rotat- 
ing objects, likely populating Box 1 or Box3. Binary evolution is 
also able to populate Box 2, as post-mass transfer systems may 
remain so ti ght that th e tidal interaction can spin down the ac- 
cretion star d Lan g er et al. 2008). 

While quantitative statements requires binary population 
synthesis models which include a methodology for non- 
conservative mass transfer (which are not yet available), it 
is clear that the canonical mass transfer evolution leads to 
a rapidly rotating, nitrogen-rich main sequence accretion star 
(iLanger et al.ll2008h . Thus, invoking binaries has the potential 
to solve the problems with Boxes 1 and 2, but for each star ap- 
pearing in these boxes, one or more stars would appear in Box 3 
due to binary interaction, as binary evolution is expected to lead 
more often into Box 3 than into Box 1 or 2. 

This implies that, using binaries to solve the problems with 
Box la and, in particular, with Box 2, worsens the problem with 
Box 3, possibly even to the extent that no rotationally mixed 
single stars are required to explain those found in Box 3. We 
conclude that close binary evolution, as far as it is currently 
understood, is unlikely to be responsible for the slowly rotating 
nitrogen-rich population of observed stars in Box 2. 

5.2. Magnetic fields 

Effects of magnetic fields other than the magnetic angular mo- 
mentum transport implemented in our models may affect a frac- 
tion of massive stars. Further down the main sequence, magnetic 
Ap stars tend to rotate more slowly than non-magnetic A stars 
and show surface abundan ce anomalies . Also the ratio of Ap 
to A stars is about 5-10% (IWolfail968b which is not far from 
the percent a ge of stars which populates Box 2. Furthermore, 
iMorel et"an (l2008h recently found magnetic fields in Galactic 
nitrogen enriched slowly-rotating early B-type dwarfs. All this 
makes the possibility that magnetic fields play a role in the pop- 
ulation of stars in Box 2 an interesting direction of investigation. 

A potential magnetic mechanism would not necessarily need 
to act simultaneously on both the nitrogen abundance and the 
stellar rotation, as is believed to occur in Ap stars. In fact one 
could envisage a mechanism that slows stars down through mag- 
netic braking, with rotational mixing occurring whilst the star is 
still a fast rotator responsible for the nitrogen enrichment. The 
fact that the number of predicted stars in Box 3 is more or less 
the sum of the stars observed in Boxes 3 and 2 may be sugges- 
tive of a migration process from Box 3 to Box 2. The lack of ob- 
served stars in Box 4 then indicates that such a migration needs 
to occur on a rapid time scale. A massive star example fo r mag- 
netic braking has been found by iTownsend et al.l (1201 Oh in the 



magnetic helium-strong star a Ori E with a characteristic spin 
down time of 1.34Myr Also, iMevnet et al.l (1201 ih have calcu- 
lated models with a simple description of magnetic braking that 
leads to slow rotation and produces a nitrogen surface enhance- 
ment by inducing strong internal differential rotation. 

On the other hand, one may also imagine a migration of stars 
from Box 5 to Box 2, e.g. through the buoyant rise of magnetic 
bubbles which carry magn etic flux and nitrogen from the con- 
vective core to the surface (iMacGregor & Cassinellill200l . We 
conclude that magnetic field effects may play a role, but they are 
not sufficiently understood to allow secure predictions. 

5.3. Comparison with earlier analyses 



5 Mvr 




Teff [kK] 

Fig. 12. HR-diagram of the 107 sample stars in the Hunter di- 
agram. Stars from the Nil and NGC 2004 fields are repre- 
sented by circles and triangles, respectively. The color coding 
represents the box into which stars fall in the Hunter diagram. 
Isochrones made from non-rotating evolution models are shown 
in green. The isochrone part below the magnitude cut is plotted 
in gray. 



Some of the conclusions deriv ed in this paper were already 
obtained by iHunter et al.l (l2008al) . based on the same observa- 
tional sample as the one scrutinized her e, but on ly based on a 
subset of the stellar models p resented in Paper J ra ther than on 
detailed population synthesis. iMaeder et al.l ( 20091) argued that 
the conclusions of Hunter et al. were unwarranted, due to the 
dispersion of the sample stars in mass, age, rotation, and bina- 
rity (we ignore the suggestion of a metallicity spread, which is 
known to be small in the LMC). In their own analysis, where 
they reduced the sample size by a factor of about 4 to consider 
only stars of similar mass, they still maintained the dispersion in 
age (see Fig. [T2l) . rotation and binarity, and found their models 
of rotating single stars to explain the observations satisfactorily. 

We find the picture emer ging from the reduced sample ana- 
lyzed bv lMaeder elaD (l2009h and that from lHunter et al.l(l2008ah 
not very dififerent. This is in fact not surprising, since the restric- 
tion in mass applied by Maeder et al. is not necessary. The mass 
range of stars in the full sample is already so small (9...25M0), 
that the expected nitrogen enrichment is very similar through- 
out the considered mass range (cf. Fig. 8, and also Fig. 1 of 
iMaeder et al.l (l2009h V This is in fact confirmed by their analysis. 
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which shows th e bands of potentiall y rotationally mixed stars in 
Figs. 3 and 4 of^Maed er et al.l (l2009l) being exactly the same, im- 
plying that the mass dependence is negligible. While Maeder et 
al. conclude that the stars in Box 3 agree with single star models, 
we do the same, only we note that binary evolution provides an 
alternative explanation for these stars. 

Maeder et al. also find stars in Box 2 in their reduced sam- 
ples, with a similar proportion as they occur in the full sample. 
While Maeder et al. argue that they may occur due to binary 
effects and thus do not disturb t he agreement of sing le star pre- 
dictions with the observations, iMevnet et al.l (1201 Ih argue that 
massive slowly rotating nitrogen-rich main sequence stars could 
be single stars unde rgoing magnetic braking, as proposed by 
iHunter et al.l (l2008al) . Only the fraction of rapidly rotating stars 
without strong nitrogen enrichment (Box 1) in the iMaeder et alJ 
(120091) sample is smaller than in the full sample, which is due 
to the small size of the samples of Maeder et al. However, of 
course, all Box 1 stars in the full sample remain to be explained. 

In any case, we developed the population synthesis approach 
as presented in this paper in order to construct a tool to ana- 
lyze also diverse stellar samples in a statistically sound way. The 
spread in age, mass and rotation (and in principle, metallicity) is 
fully accounted for by our me thod. Our new results reinforce the 
uncomfortable conclusions of lHunter et al.l (l2008al) that two sig- 
nificant groups of massive core hydrogen burning stars stand out 
as being in conflict with evolutionary models, posing a challenge 
to the concept of rotational mixing. 



6. Conclusions 

Through detailed single star population synthesis modeling, we 
have tried to reproduce the properties of the largest homoge- 
neous sample of early B-type stars in the LMC. This sample 
is the first for which nitrogen surface abundances have been 
determined for a wide range in projected rotational velocities. 
Uncertain mixing parameters in the underlying stellar evolu- 
tion models (overshooting and rotational mixing) have been cal- 
ibrated such that they reproduce the observed main-sequence 
widening deduced from t he sharp drop in rotational velocities 
(lPaDerlLIVinketal.ll2010h . as well as the observed pattern of ni- 
trogen enrichment as a functio n of the projected rota tion velocity 
in the majority of fast rotators (iHunter etal."2008a'). 

We confirm quantitatively the result of Hunter e t al.l (l2008al) 
that two sub-populations of stars, namely rapidly rotating, un- 
enriched stars and slowly rotating, nitrogen-rich stars (Boxes la 
and 2 in Fig. [TOl i. are not reproduced by a population synthesis 
model for single stars. 

While the group of rapidly rotating nitrogen-rich stars 
(Box 3) is, by construction, in qualitative agreement with our 
simulation, only about half of the number of predicted stars are 
found to be present in the observational sample. Although possi- 
bly explained by uncertainties in our method (see Sec. 14. Il l, pre- 
dictions of close binary evolution models, which also produce 
nitrogen-rich rapid rotators, tend to increase this discrepancy. 

We find it unlikely that the group of slowly rotating nitrogen- 
rich stars (Box 2 in Fig. [TOl i is produced via mass transfer and 
subsequent tidal spin-down in close binary systems, as a binary 
population which could produce these stars would likely over- 
populate the group of rapidly rotating nitrogen-rich stars (see 
Sec. 15. 11 1. At present, it appears plausible that some of the unex- 
plained slowly rotating, nitrogen-rich stars are aff'ected by mag- 
netic fields (cf. the magneti c field measurem ents for Galactic 
stars with similar properties: iMorel et al.ll2008l) . 



While a detailed physical prescription of possible magnetic 
effects in stellar models is currently not available, several binary 
effects are much better understood. Therefore, it appears promis- 
ing to perform a binary population synthesis of massive main 
sequence stars in the near future, which will give quantitative 
constraints on the relative importance of binaries in the various 
subgroups. In fact, this appears to be required in order to arrive 
at solid conclusions on the nature of the group of rapidly rota- 
tion nitrogen-rich stars, and, thereby, obtain clear evidence for 
the universality of rotational mixing in massive stars. 

Furt her light should be shed on this by the ongoing Tarantula 
Survey dEvans et al.ll2010l) . which has targeted over 1000 stars 
in the 30Doradus of the LMC. 30 Dor is the largest Hn region 
in the Local Group, containing a large number of very massive 
stars. The new survey is tailored for the detection of binarity, 
so will help sample the mass range above 30 M©, where mixing 
effects should be more significants, and therefore more easily 
detectable. 
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Appendix A: The LMC Sample 

In Fig. lA. II we show the mass and projected rotational velocity distributions of 
our sample. In the upper panels, the targets are divided into the VLT-RLAMES 
fields. The rotational velocity distributions in both fields are very similar The to- 
tal mass distribution is dominated by Nl 1 stai's at high masses, while NGC 2004 
contributes mainly to the lower-mass end. 

We divided our sample into single and binary stars, as shown in the mid- 
dle and lower panels of Fig. lA. II respectively. The color coding corresponds to 
our selection effects that are later applied to each star of our population syn- 
thesis (see discussion in Sec. 13. 3t . The single star velocity distribution is well 
represented by those stars with measured nitrogen abundances (yellow panels). 
Around 50 km/s is a peak of low gravity objects (orange), which is somewhat 
sensitive to our overshooting calibration. A value of 50km/s is also compara- 
ble to the expected the line broadening due to macro turbulence, therefore, these 
stai's could be rotating at lower velocities. 

The mass distribution of single stars with measured nitrogen abundances 
is also representative of the mass distribution of the total single star sample. 
Nitrogen abundances were not detennined for the hmited sample of the most 
massive stars owing to the problems of modeling the N iii line. The binary veloc- 
ity distribution shows three peaks (at 30, 100 and 200 km/s), which may well be 
not statistically significant given the small number of stars. 
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Fig.A.l. Di s tributio n of evolutionary masses (left) and projected rotational velocities (right) of O- and B-type stars from 
iHunter et alj ( l2009ah in the NGC 2004 and Nil VLT-FLAMES fields. Top panels: The contributions from each VLT-FLAMES 
field are shown in yellow (Nil) and blue (NGC 2004); Middle: Apparently single stars; Bottom: Radial velocity variables. The 
selection effects discussed in Sec [33] are color-coded as per the legend in each plot. For the objects in light yellow we do not have 
nitrogen measurements for various reasons. The stars shown with gold panels are those with measured nitrogen abundances that are 
the primary focus of the current analysis. Not shown in the mass histogram are two O-type stars, with masses of 62 and 83 Mq. 
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